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Recent experiments on ultracold atomic gases in an optical lattice potential have produced a Mott 
insulating state of 87 Rb atoms. This state is stable to a small applied potential gradient (an 'electric' 
field), but a resonant response was observed when the potential energy drop per lattice spacing (E), 
was close to the repulsive interaction energy (U) between two atoms in the same lattice potential 
well. We identify all states which are resonantly coupled to the Mott insulator for E ~ U via 
an infinitesimal tunneling amplitude between neighboring potential wells. The strong correlation 
between these states is described by an effective Hamiltonian for the resonant subspace. This 
Hamiltonian exhibits quantum phase transitions associated with an Ising density wave order, and 
with the appearance of superfluidity in the directions transverse to the electric field. We suggest that 
the observed resonant response is related to these transitions, and propose experiments to directly 
detect the order parameters. The generalizations to electric fields applied in different directions, 
and to a variety of lattices, should allow study of numerous other correlated quantum phases. 



I. INTRODUCTION 



Recent experiments on ultracold trapped atomic gases 
have OBened a new window onto the phases of quantum 
mattend'El. A gas of bosonic atoms has been reversibly 
tuned between superfluid and insulating ground states 
by varying the strength of a periodic potential produced 
by standing waves of laser lightfi. These experiments of- 
fer unprecedented control of the microscopic parameters, 
and allow exploration of parameter regimes not previ- 
ously available in analogous condensed matter systems. 

This paper focuses on one such "extreme" parameter 
regime. Let w be the amplitude for an atom to tunnel 
between neighboring minima of the standing laser wave, 
and U be the repulsive interaction energy between two 
atoms in the same potential well. When w is smaller 
than a value of order U, the ground state is a Mott insu- 
lator for certain values of the atomic density or chemical 
potential. In this state, the average number of atoms in 
each potential well must be an integer, nor( see Fig [l]). 
Now consider "tilting" this Mott insulatora i.e. plac- 
ing it under an external potential which decreases lin- 
early along a particular direction in space. Conceptually, 
it is useful to imagine that the atoms carry a fictitious 
'charge', and then this potential gradient corresponds to 
applying a uniform 'electric' field, E (in practice this field 
is applied by changing the position of the center of the 
atomic trapa). We measure E in units of energy, defining 
E to be the maximal drop in potential energy of an atom 
moving between nearest-neighbor minima of the periodic 
potential (the potential energy drop depends upon the 
choice of the nearest neighbor, and we choose the direc- 
tion^) along which the drop is the largest to define E). 
In almost all Mott insulators consisting of electrons or 
Cooper pairs, all reasonable electric fields that can be 
achieved in the laboratory are small enough so that the 
relation E <C w,U is well satisfied. Remarkably, in the 
new atomic systems significantly larger 'electric' fields 
are easily achievable: this paper shall discuss the regime 
E ~ U which has been explored in the recent experi- 



FIG. 1: Figures contain schematic representations of the 
Mott insulator, and of various states coupled to it. Shown 
above is the Mott insulator with no = 2. Each well represents 
a local minimum of the optical lattice potential - these we 
number 1-5 from the left. The potential gradient leads to a 
uniform decrease in the on-site energy of atom as we move 
to the right. The grey circles are the bi bosons of (1.2). The 
vertical direction represents increasing energy: the repulsive 
interaction energy between the atoms is realized by placing 
atoms vertically within each well, so that each atom displaces 
the remaining atoms upwards along the energy axis. We have 
chosen the diameter of the atoms to equal the potential en- 
ergy drop between neighboring wells — this corresponds to the 
condition U = E. Consequently, a resonant transition is one 
in which the top atom in a well moves horizontally to the top 
of a nearest-neighbor well; motions either upwards or down- 
wards are non-resonant. 



ments of Greiner et al.u. More precisely, we shall discuss 
the regime 



\U-E\,w<E,U, 



(1.1) 



while allowing the ratio (U — E)/w to be arbitrary. 

We mention, in passing, another experimental system 
whic h has been studied under conditions analogous to 
(1.1). Electron transport has been investigated in ar- 
rays of GaAs quantum dotaJ, when the voltage drop be- 
tween neighboring quantum dots (the analog of E) is at 
or above the charging energy required to make the transi- 
tion (the analog of U). However, in these systems the ex- 
cess electron energy can be dissipated away to the under- 
lying lattice, and so it appears that the threshold behav- 
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ior can be described by dissipative classical modelsi. In 
contrast, for the atomic systems of interest in the present 
paper, there is essentially no dissipation over the time 
scales of interest, and a fully quantum treatment must 
be undertaken. 

It useful to explicitly state our model Hamiltonian for 
the Mott insulator for our subsequent discussion. We 
will consider only Mott insulators of bosons, although 
the extension to fermionic Mott insulators is possible^. 
We label the minima of the periodic potential by lattice 
sites, i, and assume that all bosons occupy a single band 
of "tight-binding" orbitals centered on these sites. Let b\ 
be the creation operator for a bason on site i. We will 
study the boson Hubbard modeEfil'S 



E^e-r.n, (1.2) 



where (ij) represents pairs of nearest neighbor sites, 



rii 



(1.3) 



Yi are the spatial co-ordinates of the lattice sites (the 
lattice spacing is unity) , and e is a vector in the direction 
of the applied electric field (e is not necessarily a unit 
vector — its length is determined by the strength of the 
electric field, the lattice structure, and our definition of 
E above). We will mainly consider simple cubic lattices, 
with the e oriented along one of the lattice directions and 
of unit length. Not shown in (1.2) is an implied chemical 
potential term which is chosen so that the average density 
of atoms per site is no- We will restrict our attention to 
the case where uq is of order unity. 

Some simple key points can be made by first consider- 
ing the non-interacting case, U = 0, and also by simpli- 
fying to one spatial dimensional. For this special case, we 
can write H as 



Ho = - (wb\b e+1 + wb\ +1 b t + Eib\b t 



(1.4) 



where I is an integer labelling the lattice sites. The exact 
single-particle eigenstates of Ho can be easily obtained: 
the eigenenergies form a Wannier-Stark ladder, and the 
most important property of the wavefunctions is that 
they are all localized. Specifically, the eigenstates can 
be labeled by an integer m which runs from — oo to oo, 
the exact eigenenergies are 



e m = Em, 



(1.5) 



and the corresponding exact and normalized wavefunc- 
tions can be expressed in terms of Bessel functions: 



i> m {f) = Ji- m (2w/E); 



(1.6) 



for a derivation see e.g. Ref. [n] (their analysis is in a 
different gauge). The m'th state is localized near the site 



m, and for large \i — m\ its wavef unction decays as 

'\e-m\E 



\i>m{£) \ ~ exp 



-\i — m\ In 



(1.7) 



the decay is faster than ex pon ential, and is extremely 
rapid under the conditions ([0]). The reader should re- 
sist the temptation to imagine that a particle placed ini- 
tially at the site £ will eventually be accelerated by the 
applied electric field out to infinity. Instead, the parti- 
cle remains localized near its initial site, and undergoes 
Block oscillations with period h/E; indeed, as is clear 
from the simple form of (1.5), its wavefunction is exactly 



equal to its initial wavefunction at regular time intervals 
of h/E. The particle can escape to infinity only be a pro- 
cess of Zener tunneling to higher bands not includ ed i n 
the single band tight-binding models in (IA) and fll.2|); 
the probability of such tunneling is negligibly small in 
the experiments of interest here, and so will be ignored 
in our analysis. 

We now return our discussion to the full Hubbard 
model (1.2). As was the case in (|1.5|), the spectrum of 



this Hamiltonian is unbounded from below for E ^ 0, 
and so it does not make sense to ask for its "ground 
state" for any density of particles. Rather, guided by 
the experimental situation of Ref. ||, we are interested in 
states which are accessible from the translationally in- 
variant Mott state (with an average of no particles on 
every site) over .the experimentally relevant time scales. 
The experimento begins at E = with a Mott insula- 
tor with no particles per site, rapidly ramps up £ to a 
value of order U, and detects the change in the state. For 
w <C U, and for most values of E, the experiments dis- 
played little detectable change in the state of the system. 
We can initially understand this by a simple extension 
of the argument presented above for the non-interacting 
model Ho- Consider a 'quasiparticle' state of the Mott 
insulator, created by adding a single additional particle 
on one site, as shown in Fig |2|a. To leading order in w/U, 
the motion of this quasiparticle along the direction e is 
described by an effective Hamiltonian which is identical 
in form to Ho, but with the hopping matrix element w 
replaced w(no + 1)- So any such quasiparticle states cre- 
ated above the Mott insulator will remain localized and 
will not have the chance to extend across the system to 
create large changes in the initial state. A similar local- 
ization argument applies to the quasihole state shown in 
Fig ^fa: it experiences an electric force in the opposite 
direction, the effective hopping matrix element is wno, 
and all quasihole states are also all localized in the di- 
rection e. Indeed, it is not difficult to see that the same 
localization argument applies to all deformations of the 
Mott insulator which carry a net charge. 

The important exceptions to the above argument for 
the stability of the Mott state are deformations which 
carry no net charge. It is the primary purpose of this pa- 
per to describe the collective properties of such neutral 
states. They will be shown to yield a resonantly strong 
effect on the Mott state when E ~ U, which has been 
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FIG. 2: Notation as in Fig. [lj. (a) A quasiparticle on site 3; 
the motion of t his quasiparticle is described by the localized 
wavefunctions (1.6) but with w replaced by 3w. (b) A quasi- 
hole on site 3; the motion of this quasihole is also described 
by the localized wavefunctions (fhq) but with w replaced by 
2w. 




FIG. 3: Notation as in Fie. § (a) A dipole on sites 2 and 3; 
this state is resonantly coupled by an infinitesimal w to the 
Mott insulator in (a) when E — U. (b) Two dipoles between 
sites 2 and 3 and between 4 and 5; this state is connected via 
multiple resonant transitions to the Mott insulator for E = U. 



dramatically observed in the experiments of Greiner et 
al.a. Indeed, Greiner et al. have already identified an im- 
portant neutral deformation of the Mott state — it is the 
dipole state consisting of a quasiparticle-quasihole pair 
on nearest neighbor sites, as shown in Fig ||a. A key con- 
sequence of our discussion above is that, for w <C E (a 
condition we assume throughout), we can safely neglect 
the independent motion of the quasiparticle and of the 
quasihole along the direction of e. Only their paired mo- 
tion as dipoles will be important along e, although they 
can move independently along directions orthogonal to 
e. 

For w — 0, the dipole state in Fig ||ei differs in energy 



from the Mott state by U — E. So these states become 
degenerate at U = E, and an infinitesimal w leads to a 
resonant coupling between them. However, there are a 
large number of other states which are resonantly cou- 
pled to one of more of these states, and they also have 
to be treated on an equal footing. Indeed, it is suffi- 
cient for an given state to be resonantly coupled to any 
one other state in the manifold of resonant states for it 
to be an equal member of the resonant family — it is not 
necessary to have a direct resonant coupling to the par- 
ent Mott insulator. The reader should already notice 
that multiple dipole deformations of the Mott insulator 
(like the state in Fig |^b) are part of the resonant family. 
In dimensions greater than one, these dipole states are 
only a small fraction of the set of resonant states, as we 
will show below. We are now in a position to succinctly 
state the purpose of this paper: identify the complete set 
of states re sona ntly coupled to the Mott state under the 
conditions (1.1), obtain the effective Hamiltonian within 
the subspace of these states, and determine its spectrum 
and correlations. The results will allow us to address the 
strong response of the Mott insulatar to an electric field 
E ~ U observed by Greiner et alB, and lead to some 
definite predictions which can be tested in future exper- 
iments. 

The first step in our program is a complete descrip- 
tion of the set of resonant states. We will do this first 



for one dimension in the Section I A , and for all higher 
dimensions in Section IB. The effective Hamiltonian in 



the resonant subspace will be shown to contain strong 
correlations among its degrees of freedom, but we will 
demonstrate that these can be satisfactorily treated by 
available analytic and numerical methods in many body 
theory. Before embarking on a detailed description of 
our computation, the reader may find it useful to ex- 
amine Figs U and |] for an understanding of the origin 
of the strong correlations in the one dimensional case. 
Fig U contains only dipole states: notice that while res- 
onant dipole states can be created separately on nearest 
neighbor links, it is not possible to create two dipoles 
simultaneously on such links (as in Fig ^a) without vi- 
olating the resonant conditions. This implies an infinite 
repulsive interaction between nearest neighbor dipoles in 
the effective Hamiltonian. Two (or more) dipoles can be 
safely created when they are further apart, as shown in 
Fig ||b. Thus the dipole resonances are not independent 
of each other, and the wavefunction contains non-trivial 
'entanglements' between them. 



A. One dimension 

It is not difficult to see that, in one spatial dimen- 
sion, the set of all nearest-neighbor dipole states consti- 
tute the entire family of states resonantly coupled to the 
Mott insulator in Fig [l] for U = E and an infinitesimal 
w. The only subtlety concerns states like those in Fig [|b, 
which are not made up of nearest-neighbor dipoles. For 
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FIG. 4: Notation as in Fig. |l|. Two states which are not part 
of the resonant manifold, (a) An attempt to create dipoles 
between sites 2 and 3 and also between sites 3 and 4; the result 
is a single dipole of length 2 which has energy U — 2E relative 
to the Mott insulator, and so this long dipole is not part of the 
resonant family of states. (6) A state with energy 3(U — E) 
relative to the Mott insulator; this state is not part of the 
resonant family because its largest effective matrix element to 
any stat e in the resonant family is of order w 2 /U (for U = E; 



1.8)). In contrast, all states within the resonant family 



are connected to at least one other state also in the family by 
a matrix element of order w. 



w = 0, this state has energy 3(C7 — E) relative to that 
in Fig [|. However reaching the state in Fig ||b from any 
state in the resonant family requires a detour through a 
non-resonant state. A simple second-order perturbation 
theory calculation shows that the closest state from the 
resonant family connected to Fig is a state with dipole 
between sites 3 and 4, and that the effective matrix ele- 
ment between them is 



w 2 n y/n (n + 1) 



1 



2U- E 



(1.8) 



this is negligibly small, under the conditions (1.1), com- 
pared to the non-zero matrix elements (= w) between 
states within the resonant family. Hence we can safely 
neglect the state in Fig [|b. More completely, the argu- 
ment is that after we diagonalize the Hamiltonian within 
the resonant family, states coupled to that in Fig ||b will 
diffe r from it by an energy of order u>] the coupling in 
(1.8) will then be too weak to induce a resonance. 

It is convenient now to introduce bosonic dipole cre- 
ation operators, d\, to allow us to specify the resonant 
subspace and its effective Hamiltonian. Let |Mno) be the 
Mott insulator with no particles on every site (the state 
in Fig |] is |M2)). We identify this state with the dipole 
vacuum |0). Then the single dipole state is 



dj|0> = 



\/n {n + 1) 



hbl +1 \Mn ) 



(1.9) 



Notice that we have placed the dipole operator on the 
left edge of the dipole which actually resides on links 
between the lattice sites. Clearly, we cannot create more 
than one dipole resonantly on the same link: hence the 
dipoles satisfy an on-site hard core constraint 



did, < 1. 



(1.10) 



Moreover, we cannot create two dipoles simultaneously 
on nearest neighbor links — this leads to a non-resonant 
state like that in Fig ||a; such states are prohibited by a 
hard core repulsion between nearest neighbor sites 



d\dtd\ +1 dt + i 



0. 



(1.11) 



The resonant family of states can now be completely 
specifie d as the se t of all states of the boson di which 
satisfy ( |1.10j) and (1.11). A typical state is sketched be- 
low in Fig [5|a. Notice that the dipole vacuum, |Mnn) is 
one of the allowed states. 

It is now a simple matter to write down the effective 
Hamiltonian, Tid for the d,. It costs energy U — E to 
create each dipole, and each dipole can be created or an- 
nihilated with an amplitude of order w (this corresponds 
to the horizontal motion of particles in Figs |l|-^|) . So we 
have 

H d = -uVM"o + 1) { dl + d l) +( U - E )J2 d ^ de - 

i i 

(1.12) 



The Ham iltonian ( 1.12] ), along with the constraints 



( 1.10 , 1.11 ), constitute one of the correlated many-body 
problems we shall analyze in this paper. The eigenstates 
of Tid are characterized by no and the single dimension- 
less number 



A 



U-E 



(1.13) 



and a description of their properties as A ranges over all 
real values is in Section ||. (Strictly speaking, the eigen- 
states of Hd depend only X/y/ no (no + 1), but A and no 
do not combine into a single constant in higher dimen- 
sions.) 

It is interesting to note that there is no explicit hopping 
term for the di bosons in Tid'. it appears that the bosons 
only only allowed to be created from, and t o disappear 
into, the vacuum by the first term in ( 1.12|) . However, 
this is misleading: as we will see in Section |H| the combi- 
nation of the terms in (1.12) and the constraint (1.11) 
does gen erate a local hopping term for the di bosons 
(see (2.1)). Additional dipole hopping terms also arise 
from virtual processes of order w 2 /U in the underlying 
Hubbard model 7i; however, these are negligibly small 
compared to those just mentioned and do not need to be 
included in Hd- 

We close this subsection by noting that the Hamilto- 
nian Tid in ( 1.12p and the constraints (1.10,1.11) can also 
be written in the form of a quantum spin chain. We iden- 
tify the dipole present /absent configuration on a site £ as 



a pseudospin cr| up/down (<j x ' v ' z are the Pauli matrices). 
Then cr| = 2dldg — 1 and 

H d = [-WM^o + IK + {U- E)(af + l)/2 



+ J(<7f + l)(<7f +1 + l) 



(1.14) 



The constraint (1.11) is implemented by taking the J — > 
oo limit of the last term. The spin chain model so ob- 
tained is an S — 1/2 Ising spin chain in both transverse 
and longitudinal fields. This is known not to be inte- 
grable for finite J, but it does appear that the problem 
simplifies in the J — > oo limit we consider here. 



B. Higher dimensions 

We consider here only hypercubic lattices in D spatial 
dimensions, with e oriented along one of the principal 
cubic axes and a lattice spacing in length (e.g. D = 3 
and e = (1,0,0)). Other lattices, and other directions of 
e, also allow for interesting correlated phases and these 
will be mentioned in Section IV. 



Extension of our reasoning above quickly shows that 
the dipole states now constitute only a negligibly small 
fraction of the set of states in the resonant family. Once a 
dipole has been created on a pair of sites separated by the 
vector e, its quasiparticle and quasihole constituents can 
move freely and resonantly, with matrix elements of or- 
der w, in the (D — 1) directions orthogonal to e. Allowing 
this process to occur repeatedly (while maintaining some 
constraints discussed below), we can build up the set of 
all resonantly coupled states. A typical res ona nt state in 
D = 2 is shown in Fig [s]b. As in Section I A, it is use- 
ful to give an operator definition of the resonant family. 
To allow us to distinguish between the directions paral- 
lel and orthogonal to e, we replace the £)-dimensional 
site label i, by the composite label (£,n), where I is an 
integer measuring the co-ordinate along e (as in the one- 
dimensional case), while n is a label for sites along the 
(D — 1) transverse directions. Rather than using dipole 
operators, we now want to work with bosonic quasipar- 
ticle (jp\ n ) and quasihole (h\ n ) operators, which create 
states like those in Fig and Fig ||b respectively. More 
precisely, we now identify |Mno) with quasiparticle and 
quasihole vacuum |0), and so 
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IMnn 



h\ n \Q) = -L&f,„|Mn >. 



(1.15) 



The set of resonant states can now be specified by a 
few simple constrain ts on these operators, which are the 
analogs of (1.10,1.11). First, there are the obvious on-site 
hard-core constraints that no site can have more than one 
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FIG. 5: Typical states in the resonant subspace for (a) D = 1 
and (b) D — 2. Black circles represent sites with quasipar- 
ticles (these sites have m = no + 1 (see ([L.3|)), grey circles 
represent quasiholes (these sites have nj = no — 1 ), wh ile 
the remaining sites have rii = no. Note that Qe in (1.17) is 
zero for each column i. e. the total number of quasiparticles in 
every column equals the total number of quasiholes in the col- 
umn to its immediate left. Only in D = 1 does this constraint 
imply that all states contain only nearest-neighbor dipoles. 



particle or hole: 



Pi tV pt, n < 1 
h\ n ht, n < 1 



(1.16) 



Additionally, because of the manner in which these quasi- 
particles and quasiholes appear from the Mott state, the 
total number of quasiparticles in the D — 1 dimensional 
layer with co-ordinate t + 1 must equal the total number 
of quasiholes in layer £: 

Qe = y^ (pi+i,nPe+i,n - h\ n ht n ^ = 0. (1.17) 

n 

While the quasiparticles and quasiholes are allowed to 
move freely within each D— 1 dimensional layer, they can- 
not move resonantly out of any layer on their own; this is, 
of course, related to the localization of the Wannier-Stark 
ladder states discussed earlier in this section. 



Continuing the analogy with Section [A, we can now 
easily write down the effective Hamiltonian, H. p h, for the 
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quasiparticles an d qua sihol es wh ich acts on the set of 
states defined by ( 1.1 6| ) and ( 1.17 ). The terms in th e firs t 
two lines are the same as those already present in (1.12), 
but expressed now in terms of the quasiparticle/hole op- 
erators, while the last line is associated with motion along 
the transverse D — 1 directions: 

T-tph = -wy/n (n + 1) y ] (p£+i t nhi t „ + p\ +lM h\ M ) 

^^E(^n + 4n^n) (1-18) 



+ 



w 



^ (n h\ n hi, m + (n + ^)p\ n Pt, m + h.c.) . 

£, {nm) 

Here (nm) represents a nearest neighbor pair of sites 
within a single (D — 1) dimensional layer orthogonal to e. 
Notice that all the Qi in (1.17) commute with TL p h, as is 
required for the consistency of our approach. As was the 
case in one dimension, the properties of H p h are deter- 



mined by the single dimensionless constant A in (1.13); 
these will be described in Section |TTT[ 

It is worth reiterating explicitly here that upon special- 
ization to the case of D = 1 (when the indices n, m only 
have a single allowed value and the set (nm) is empty), 
the Hamiltonian TL p h above is exactl y equ ivalent to the 
one dimensional dipolc model Hd in (1.12). 

We note in passing that in a manner similar to TLd, 
H.ph in (1.18) can also be written as a S = 1 spin model, 
with the empty/qausiparticlc/quasihole states on a site 
corresponding to spin states with S z = 0,1, —1. 

The outline of the remainder of the paper is as follows. 
The properties the D = 1 model Tld will be described in 
Section O, while the D > 1 model TL p h will be consid- 
ered in Section III . We discuss extensions of our results to 



other lattices and field directions in Section IV. Implica- 
tions of our results for experiments appear in Section [v|. 
The appendices contain some technical discussion on the 
nature of the quantum phase transitions found in the 
body of the paper. 



II. DIPOLE PHASES IN ONE DIMENSION 

This section will describe the spect rum of the one- 
dimensional dip ole H amilto nian TLd hi (1-12), subject to 
the constraints ( 1.10 ) and ( 1.11 ). 

An essential point becomes clear simply by looking at 
the limiting cas es A — > oo and A — > — oo (the coupling A 
was defined in ( 1.13| )). For A — » oo the ground state of 
TLd is the non-degenerate dipole vacuum |0). In contrast, 
for A — * — oo the ground state is doubly degenerate, be- 
cause there are two distinct states with maximal dipole 
number: (• • • d\d\d\ ■ ■ • )|0) and (• • • 444 ' " ' )l°>- This 
immediately suggests the existence of an Ising quantum 
critical point at some intermediate value of A, associ- 
ated with an order parameter which is a density wave of 
dipoles of period two lattice spacings. We will shortly 



present numerical evidence which strongly supports this 
conclusion. 

Further analytic evidence for an Ising quantum crit- 
ical point can be obtained by examining the excitation 
spectra for the limiting A regimes, and noting their sim- 
ilarity to those on either side of the critical point in the 
quantum Ising chainO. 

For A — > oo, the lowest excited states are the single- 
dipoles: \£) = dt\0); there are TV such states (N is the 
number of sites), and, at A = oo, they are all degenerate 
at energy U — E. The degeneracy is lifted at second order 
in a perturbation theory in 1/A: by a standard approach 
using canonical transformations, these corrections can be 
described by an effective Hamiltonian, TLd, c s, that acts 
entirely within the subspace of single dipole states. We 
find 



H cl , eS = (U-E)J2[\t)(t\ 



n (n + 1) 
A 2 



(\i)(£\ + \i)(e + i\ + \i+i)(e\) (2.i) 



Notice that, quite remarkably, a local dipole hopping 
term ha s ap peared, as we pro mised ear lier at the end of 
Section [I A[ The constrain ts ( J l.10 1.11 ) played a crucial 
role in the derivation of ( |2.l| )l Upon consi dering per- 
turbations to \£) from the first term in (1.12) it initially 
seems possible to obtain an effective matrix element be- 
tween any two states \£) and \£'}. However this connec- 
tion can generally happen via two possible intermedi- 
ate states, \£) -> d\d\, |0) -> and \£) -> |0) -> \i'), 
and the contributions of the two processes exactly can- 
cel each other for most £, £'. Only when the constraints 



( 1.10 , 1.11 ) block the first of these processes is a residual 
matrix element possible, and these are shown in ([2.l|). 
It is a simple matter to diagonalize TLd, e s by going to 
momentum space, and we find a single band of dipole 
states. The lowest energy dipole state has momentum 
7r: the softening of this state upon reducing A is then 
consistent with the appearance of density wave order of 
period 2. The higher excited states at large A consist of 
multiparticle continua. of this band of dipole states, just 
as in the Ising chainO. 

A related analysis can be carried out for A — > — oo, 
and the results are very similar tp-,those for the ordered 
state in the quantum Ising chairJlil. The lowest excited 
states are single band of domain walls between the two 
filled dipole states, and above them are the corresponding 
multiparticle continua. 



A. Exact diagonalization 

We numerically determined the exact spectrum of TLd 
for lattice sizes up to N — 18. As will be evident below, 
these sizes were adequate to reliably extract the limiting 
behavior of the N — > oo limit. 

The complete spectrum of Tld is shown in Fig ^| for 
N = 8 and no = 1. We used periodic boundary con- 
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FIG. 6: All the eigenvalues of TLd for TV = 8 and no = 1. 
Note that the ground state is non-degenerate for positive A, 
and there are two low-lying levels with an exponentially small 
splitting for A < and |A| large. 



ditions on the dipolc Hamiltonian in (1.12). Note that 
these do not correspond to periodic bound ary conditions 
for the original model (1.2); indeed, for (1.2) the pres- 
ence of the electric field implies that periodic boundary 
conditions are not physically meaningful. Nevertheless, 
it is useful to apply periodic boundary c ondit ions to the 
translationally invariant effective model ( 1.12| ) , merely as 
a mathematical tool for rapidly approaching the N — ► oo 
limit. Note that Fig ^ shows a unique ground state for 
A — * oo and a two- fold degenerate state for A — > — oo. 
Above these lowest energy states, there is a finite energy 
gap, and the excited states have clearly split into bands 
corresponding to the various "particle" continua; these 
"particles" are dipoles for A — > oo, and domain- walls be- 
tween the two ground states for A — > — oo, as we discussed 
in the perturbative analysis above. 

We test for a quantum critical point at intermediate 
values of A by plotting the energy gap, A, in Fig 0. This 
gap is the spacing the between the lowest two of the 
eigenvalues plotted in Fig || (for finite system sizes, these 
low- lying levels are always non-degenerate). It becomes 
exponentially small in the system size as we approach 
the two degenerate ground states which are present for A 
sufficiently negative. In the opposite limit, A approaches 
a finite non-zero value, which becomes U — E, for A large 
and positive. If these two phases are separated by a quan- 
tum critical point, we expect the energy gap to scale as 
A ~ N~ z at the critical point A = A c , where z is the 
dynamic critical exponent. The Ising critical point has 
z = 1, and so Fig |^ plots iVA as a function of A. We 
observe a clear crossing point at A c —1.850 which we 
identify as the position of the Ising quantum phase tran- 
sition. Note that the critical point is shifted away from 
the naive value E — U (A = 0) to E > U because of 
quantum fluctuations associated with the hopping of the 
dipoles. 
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FIG. 7: The spacing between the lowest two eigenvalues of Tid 
(= A) as a function A for various system sizes and no = 1. 
We used periodic boundary conditions for Tid- 
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FIG. 8: Scaling plot of the energy gap to test for (2.2). We 
used A c = —1.850 and no = 1. 



A second test of Ising criticality is provided by also 
rescaling the horizontal axis of Fig ^| with N. General 
finite size scaling arguments imply that the energy gap 
should obey the scaling form 

A = /V-^(/V 1/!/ (A-A c )) (2.2) 

where is a universal scaling function, and v is the cor- 
relation length exponent. We test for (2.2) in Fig 0with 



the Ising exponent v = 1, and again find excellent agree- 
ment. 
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FIG. 9: Scaling plot of numerical res ults for the order param- 
eter structure factor, SV, defined in (2.3). We used no = 1. 



A final, and most sensitive, test for Ising criticality is 
provided by a measurement of the anomalous dimension 
of the order parameter. The order parameter is the den- 
sity of dipoles at momentum n, and so we computed its 
equal-time structure factor 



(2.3) 



Standard scaling arguments imply that this should scale 
as N 2 ~ z ~ v at A = A c , where r\ is the anomalous di- 
mension of the order parameter. Using the Ising expo- 
nent rj — 1/4, we expect SV ~ TV 3 / 4 . This is tested 
in Fig |H Note that there is an excellent crossing point 
at A c w —1.853. This position of the crossing point is 
completely consistent with the crossing point found in 
Fig 0. Thus Fig ^| provides strong evidence for the ex- 
pected Ising exponent r\— 1/4. We have also examined a 
plot which scales the horizontal axis in Fig || as in Fig [|: 
the data collapse is again excellent. 



III. QUASIPARTICLE AND QUASIHOLE 
PHASES IN HIGHER DIMENSIONS 

This section will discuss the properties of the D > 1 di- 
mensional model of the p^ n quasiparticl es and h( n quasi- 
holcs described the Hamiltonian H p h in (LIS), subject to 
the constraints (|1_.16 1.17). 

As in Section |n it is instructive to first look at the two 
distinct limiting values of A. The nature of the ground 
states is very similar to those in D = 1 for these ranges 
of A. For A — > oo, we have a unique ground state which 
contains only small perturbations from the quasiparticle 
and quasihole vacuum |0). For A —> — oo, it is clear that 
we want to maximize the total number of quasiparticles 
and quas iholcs in the ground state, subject to the con- 
straints (1.16 1.17). There are a very large number of 



ways of doing this, but by considering perturbative cor- 
rections to the ground state energy in powers of 1/|A|, 



it is not difficult to see that it pays to choose one of 
two regular arrangements, in which the occupation num- 
bers are independent of n: ]J n ]J e cvon p\ +1<n h\ <n \0) or 

LLIL odd^+i.rXrJ )- So there is a two-fold degener- 
ate ground state for A < and |A| large, associated with 
a broken translational symmetry and the development of 
density wave order of period 2 in the longitudinal direc- 
tion, in both the quasiparticle and quasihole densities. 

The excitation spectrum in the limiting ranges of A 
can also be determined as in Section [n|. However, the 
computations are more involved and we limit ours elves 
to an analysis of the A — > oo case in Section [II A . We 



will investigate physics at intermediate values of A in the 
subsequent subsections, where we will see that the pos- 
sibilities are richer than the appearance of a single Ising 
quantu m crit ical point between the states just discussed: 
Section IIIB will present a mean field theory, while Ap- 
pendices |A| and ^ will discuss continuum quantum field 
theories which can describe long-wavelength fluctuations 
near the phase boundaries. 



A. Excitations for A large and positive 

There is a large manifold of lowest excited states, all 
of which have energy U — E, in the limit A — ► oo. These 
are the states with exactly one p quasiparticle and one 
h quasihole, with the particle on the D — 1 dimensional 
layer £ + 1 and the hole on the layer I. We label these 
states by 



\£,n,m)=p\ +l h\ |0) 



(3.1) 



We break the degeneracy between these states by consid- 
ering corrections in pow ers of 1/ A. At order 1/ A, the term 
in the last line in ( 1.18 ) will allow the quasiparticle and 
the quasihole to hop independently in their own layers, 
but will not induce any couplings between states with 
different values of I. The latter appear at order 1/A 2 , 



when as in (2.1), a nearest neighbor dipole pair can hop 
longitudinally between neighboring layers; again, as in 
D = 1, the constraints (1.16 1.17) play a crucial role in 
determining these perturbative corrections. These pro- 
cesses are described the following effective Hamiltonian 
for the manifold of excited states with energy ss (U — E): 



n ph , eS = (u-e)Y, 



^2 %n,m) (£,n,'i 



(nm) ,k 



(n + l)\£,k,n){£,k,m\ 



n (n + 1) 
A 2 



^2(\£,n,ri) (£,n,', 



+\£,n,n)(£+ l,n,n\ + \£+ l,n,n)(£,n,n\ 



(3.2) 
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Note that the first summation over n, m is unrestricted 
and ranges independently over the two variables, while 
the second is over nearest neighbor pairs (nm) . 

The Hamiltonian 7i p / li0 ff can be analyzed by the stan- 
dard techniques of scatterin g th eory. The terms within 
the first two summations in (gj) lead to a "two-particle" 
continuum of quasiparticle and quasihole states, while 
the terms within the last summation allow these states 
to scatter and possibly form a dipole bound state. We 
first form states with total transverse momentum Q^, 
and relative transverse momentum qj_ (these momenta 
are D — 1 dimensional vectors) 

n,n 

(3.3) 

where Nj_ is the number of sites in each layer, and r n are 
the spatial positions of the sites. In this basis of states 
Ti-ph,e« Next, we also transform the single longitudinal 
co-ordinate, i, to a 'dipole momentum', gy: 



kn,Q_L,q_L> = 



> = -E 

Nil ^ 

ii e 



e^^Q^qj 



(3.4) 



In this basis of states, H. p h, e g takes a form which makes 
the mapping to standard scattering theory very explicit. 
The total transverse momentum, Q^, and the longitudi- 
nal dipole momentum, q» , are conserved, while there is 
scattering between different values of : 



[£ P (q±) +Sh(Q± - q±)] 



x |?||,Q±,q±) (<7||,Q±,q±| 
w 2 n (n Q + 1)(1 + 2cosg||) 
N ± (U-E) 

x kn>Qj->qj-><9ii)Q-Ljqx| > ( 3 - 5 ) 



where (for a hypercubic lattice) 
(U-E) 



£ P (qj-) 



2w(n + 1) cos(q± a ) (3.6) 



(3_6 

(3 



and the summation over a extends over the D — 1 com- 
ponents of q^. The expression for Eh((l±) is identical to 
but with uq + 1 replaced by no. The Hamiltonian in 
is that of a particle moving in (D — 1) dimensions 
with momentum q± and dispersion e p (q±) + Sh(Qj_ — 
qjj, scattering off a delta function potential at the ori- 
gin with strength w 2 no(no + 1)(1 + 2 cosq\\)/(U — E). 
Its solution is well known: in addition to the scattering 
states, a bound state must be present in D — 1 = 1,2 
for any infinitesimal attractive potential, and for strong 
enough attraction for D — 1 > 2. So for the physically 
relevant cases of D = 2,3, a bound state must form for a 
range of g|| values near ir. It is clear that the lowest en- 
ergy bound state has Qj^ = and g|| = tt: this is a dipole 



state, and just as in D = 1, it is a precursor to the ap- 
pearance of longitudinal density wave order of period 2. 
The appearance of this dipole bound state suggests that 
the first quantum phase transition out of the featureless 
and gapped phase present for large positive A is into a 
state with Ising charge order; however, our discussion 
here is for a system with a well-developed gap to quasi- 
particle and quasihole states, it is yet not clear whether 
this approach continues to hold when the gap becomes 
small — we will return to this question in Appendix |A[ 



B. Mean field theory 

This section will present the results of a mean field 
analysis of Ti. p h ■ The central idea of the mean field the- 
ory is very simple: we treat the quantum fluctuations 
along the longitudinal direction for all n by the exact 
numerical treatment developed in Section [ij for D = 1, 
while the transverse couplings are treated in a mean field 
manner. One important b enefit of this approach is that 
the important constraints (1.16) are treated exactly. 

This approach also naturally suggests the appearance 
of additional phases which have no analog in the D = 1 
case. In particular, the motion of single p and h bosons 
in the transverse direction implies that superfluid order 
can develop along these D — 1 dimensions only. There is 
no possibility of superfluidity in the longitudinal direc- 
tion because motion along this direction can occur only 
via char ge neu tral dipole pairs which appear in the first 
term in ( 1.18| ). Thi&-Lransverse superfluid therefore has 
a 'smectic' characteiO, and its existence implies that we 
have to allow for (p) and (h) condensates: these appear 
naturally in our mean field theory. 

As in the meaft,field treatment of the zero field boson 
Hubbard modelaQ, the approximation involves a decou- 
pling of a hopping term. In particular, we only decouple 
the last transverse hopping term in (1.18), and obtain 
the following mean field Hamiltonian for a set of sites, 
labelled by £, representing any chain along the longitudi- 
nal direction 



H P h,mf [<Pe>, (he)} = 
— wZ 



wZtiq (^(hi)h\ + (fi£ 

Tlo + 1) ({Pi)p\ + {Pl)*Ptj 

wy/n (n + 1) (pe+ihf + p\ +1 h\^ 
^ ~ — (p\pt + h\h^j 



( 



(3.7) 



Here Z is the co-ordination number of any site along the 
D — 1 transverse directions, and the expectation values 
(he) and (pe) have to determined self-consistently from a 
diagonalization of (3.7) subject to the constraints associ- 
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ated with (1.16), which now become 

pipe < 1 ; h\ht < 1 ; p\pih\h t = 0. (3.8) 



We have imposed the constraints ( pL . 1 7|) in a mean-field 
manner by chemical potentials /i^, whose values must be 
chosen so that 



(pl+iPf+i) = (h\he) 



(3.9) 



is obeyed; note that these constraints are macroscopic, 
and so there is no approximation involved in using a 
chemical potential to impose them. In practice, the di- 
agonalization of H p h,mf [ipe)j (he)] must be carried out 
for a finite number of sites I = 1 . . . N; we found that 
the mean-field solutions approached the N = oo limit at 
quite small and manageable values of N. The ground 
state energy of H. p h per chain is not equal to the lowest 
eigenvalue, E , of H p h, m f [(Pe), (he)] but is easily com- 
putable from it by the relation 



Eph.rnf = E + WZJ2 M(M| 2 + ("0 + l)|<P<}f 

I 



(3.10) 



Indeed, we can also view (pe) and (he) as arbitrary com- 
ple x num bers which are determined by the minimization 
of J3J0|) . 

The results of the solution of the above mean-field 
equations are summarized in the schematic phase dia- 
grams in Fig |Io| and in the numerical results in Figs [TT| 
and [l2| It is useful to discuss the phases, in turn, as a 
function of decreasing A. 

For A very large and positive (to the right of the point 
A in Fig 10), no symmetry is broken, and we have a fea- 
tureless ground state with no superfluidity and an energy 
gap to all excitations. 

There is a superfluid-insulator transition at A driven 
by the condensation of the p and h bosons. The super- 
fluidity appears in the direction transverse to the electric 
field, and all layers behave equivalently. We will examine 
fluctuations near this critical point in Appendix ^ and 
show that the interlayer coupling is irrelevant near the 
critical point in D — 3, and so each layer is described by 
an independent critical theory. 

As shown in Fig [To|a, the mean-field theory exhibits 
a first order quantum transition at the the point C as- 
sociated with the sudden development of Ising density 
wave order i.e. the states with A to the left of C have 
(PevenPeven) + (P G dd p odd)' and similar for the density 
of the h bosons. In mean-field theory, the state to the 
immediate left of C also has the loss of the p conden- 
sate in the odd layers (say), and the loss of h condensate 
in the even layers. In general, it is quite possible that 
fluctuations, beyond those included in the present mean- 
field theory, will replace the first order transition at C by 
two second order transitions at C\ and C%, as shown in 
Fig |Io|b. At the first transition at Ci, the order param- 
eter is only the Ising density wave, while there is p and 
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FIG. 10: Schematic phase diagrams as a function of A. In 
(a) we display the topology of the phase diagram found by 
the solution of the mean field equations: thin lines are sec- 
ond order quantum phase transitions, while the thick line is 
a first order transition. The parity of the I index is indicated 
as a subscript to the p operators. The expectation values of 
h quasiholes obey the same relations as those for the p quasi- 
particles, but with the roles of 'even' and 'odd' interchanged. 
The Ising density wave order is present for all A to the left of 
C . In (6) we display a hypothetical phase diagram, possibly 
induced by fluctuations, in which the first order transition is 
replaced by two second order transitions; now Ising order is 
present at A to the left of Ci. There are superfluid-insulator 
transitions at A, B, C and C2, and Ising density wave tran- 
sitions at C and C\. 



V 0.1 




15.0 



FIG. 11: Mean field numerical values of the condensates (pe) 
and (he) as a function of A for no = 1 and Z = 4. The 
solutions shown are obtained by diagonalizing (3.7) for N = 6, 
but essentially identical results obtain for N — 4. 
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FIG. 12: As in Fig [TTJ. The values of (h odd + fl even ) are very 
close, but not identical, to the p values shown above. 



h transverse superfluidity in all the layers; the second 
transition at Ci involves the continuous vanishing of the 
p (h) condensate in the odd (even) layers in a supcrfluid- 
insulator transition, in the presence of a background of 
Ising density wave order. 

The final transition at the point B involves loss of all 
p and h condensates. There is long-range Ising density 
wave order at all A to the left of B, and a gap to all 
excitations. 

The theory of fluctuations about these mean-field re- 
sults is discussed in Appendices |A| and [b|. As we have al- 
ready noted, these could be strong enough to also modify 
the topology of the phase diagram in Fig |lO|a. One ex- 
treme possibility is that the transverse superfluid phases 
could disappear entirely, and we are left only with two 
insulating phases, one with Ising density wave order and 
the other without; the phase diagram is then as in D = 1. 
However, we show in Appendix H that the value of a par- 
ticular critical exponent determines that this is not the 
generic situation. 



IV. OTHER FIELD ORIENTATIONS AND 
LATTICES 

Our discussion has so far limited itself to hypercubic 
lattices, with the direction of the electric field, e, oriented 
along one of the principal axes. Similar analyses can be 
carried out for other lattices and for other directions of 
e. A large variety of correlated phases appear possible, 
including many not related to those already discussed. 
We will illustrate these possibilities by an example here, 
but leave a more detailed discussion to future work. 

Consider a square lattice (in D = 2) but with e = 
(1,1). In this case, the resonant transitions from the 
Mott insulator involve moving a bi boson by one lattice 




FIG. 13: A typical state in the resonant subspace for a square 
lattice with e — (1, 1). Representation is as in Fig|. The 
quasiparticles and quasiholes occur only occur in dipoles ori- 
ented along the +x or +y directions. Note that it is possible 
for dipoles to undergo a ring-exchange around a plaquette, in 
which the configuration around plaquette a can become like 
that around plaquette b; this pro cess is contained in the reso- 
nant model Ti'd in ( [f.l| ) and (4.2), and does not require much 
weaker virtual processes in the Hubbard model Ti (which are 
suppressed by powers of w/U). 



spacing, either along the +x or the +y direction. How- 
ever, once such a dipole has been created, the quasipar- 
ticle and the quasihole cannot move resonantly to any 
other sites (except by processes of order w 2 /U which 
we have consistently neglected here). So the resonant 
subspace can be described completely in terms of dipole 
states, just as in the D — 1 case discussed earlier. A 
typical state is illustrated in Fig O. The effective Hamil- 
tonian of this space of dipole resonant states is identical 
in form to ( 1.12| ): 



n' d = -Wv^K+T)^ (d a + 4) + (u-e)J2 4d a , 

a a 

(4.1) 

except now the label a extends over the links of the 
square lattice. There cont inues to be a hard-core con- 
straint d\d a < 1 like ( l.lOj ), but the possibility for new 
physics arises fro m the complexity of the generalization 
of the constraint (1.11), which is now 



d a d a d\dh — for links a, b which share a common site. 

(4.2) 

Note that each dipole blocks the occupancy of dipoles 
on six neighboring links. It would be interesting to de- 
termine the properties of Tl' d subject to the constraint 
©■ 

The possibility of rich physics becomes apparent in 
thinking about the case A < and |A| large. Here the 
low-lying manifold of states corresponds to maximizing 
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the number of dipoles, and these are in one-to-one corre- 
spondence with the close-packed dimer coverings of the 
square lattice. A natural ring-exchange term of the dipole 
bosons also becomes apparent upon considering pertur- 
bative corrections in powers of 1 / 1 Ah this derivation is 
similar in spirit to that in Section y (see Fig [t3]). We 
emphasize that the dominant ring exchange does not 
come from virtual higher order processes in the under- 
lying Hubbard model Ji (which are strongly suppressed 
by factors of w/U), but is already contained within the 
physics of the resonant subspace as described by ( }4.1| ) and 
(4.2). In. analogy with other studies of quantum dimer 
modelscJO'Eiil and boson ring exchange models, possibili- 
ties of bond-ordered,phases open up. Fractionalized, and 
Bose metal phascsEj are also possible, but these may be 
more likely on non-bipartite lattices. 

We close by noting that it is easily possible to orient 
e so that only one direction is resonant. For a cubic lat- 
tice in D = 3 this can be done by choosing e = (1, a, b) 
where a, b ^ 0, 1 are some arbitrary real numbers. Then 
resonant transitions to dipole states can occur only along 
the x direction, and the resonant manifold separates into 
decoupled one dimensional systems, each of which is sep- 
arately descri bed b y the one-dimensional dipole Hamil- 
tonian TLd in (1.12). This may be a simple way of exper- 
imentally realizing the model TLd- 



V. IMPLICATIONS FOR EXPERIMENTS 

An important issue that must be faced at the outset is 
the extent to which the non-equilibrium time-dependent 
experiments can be described by the ground and low en- 
ergy states of the effective models that have been dis- 
cussed in this paper. In the experiments of Greiner et 
al.a, the 'electric field' (in practice, this is realized by 
a magnetic field gradient) is turned from an initial zero 
value to E in a time of order h/w. In a system under 
the conditions (1.1), this may not allow easy access of 
the ground state. As an alternative, we suggest that E 
be ramped up rapidly to a value to the right of the point 
A in Fig [h], and then slowly increased through the pos- 
sible critical points in Fig [K]. This could produce states 
with either the density wave Ising order, or the transverse 
superfluid order. 

Having produced such states, the next challenge is to 
directly detect the quantum order parameters associated 
with the phases in Fig [l0| We address two possible 
probes in the subsections below. 



A. Momentum Distribution 

One experimental quantity that is relatively easy to 
measure is the momentum distribution of the atoms con- 
tained in the optical lattice. This is done by shutting off 
the lattice potential and the trapping potential and allow- 
ing the atoms to freely expand until the resulting cloud 



is large enough that its density profile can be spatially 
resolved optically. The scale to which the cloud expands 
before measurement can be made much larger than the 
original lattice dimensions. In this limit the final spatial 
position at which an atom is detected determines the mo- 
mentum at which the momentum distribution function is 
being measured. 

The momentum distribution for the boson Hubbard 
model containing N sites is given by 



n(q) = |/(q)| 2 ^ e lq - ( ^- Ffc 



b% 



(5.1) 



where /(q) is the form factor for the tight- binding or- 
bitals associated with the lattice potential, and the mo- 
mentum q = rriR/(M ex ), where R is the distance from 
the detection position to the center of the trap, m is the 
mass of the atoms, and t ex is the time elapsed in the 
expansion (this expression ignores the influence of grav- 
ity, but an appropriate modification is straightforward). 
The development of off-diagonal long-range order peaks 
the momentum distribution at the values of q equal to 
the reciprocal lattice vectors of the optical lattice poten- 
tial, and has been used as an experimental signature of 
the superfluid phase. Etc 

Let us first consider the D — 1 case. A very important 
consequence of our restriction to the subspace of resonant 

states is that the boson correlator (btbp) vanishes for 



\£— 1'\ > 1. Hence Eq. (5.1) becomes (q is the component 
of q in the direction of the 'electric' field 



n 1D (g) = |/(g)| 2 



n 



^n (n + 1) 



,(5.2) 



where the lattice spacing has been taken to be unity, and 



d\ is the dipole creation operator defined in (1.9). For 



the periodic boundary conditions we have used (as we 
noted earlier, such boundary conditions are not physi- 
cal, but they should not modify the results in the limit 
of large system sizes), the values of (dg) depend only on 
the parity of i (a very small ordering field is applied to 
lift the Ising symmetry, and choose one of the ground 
states in the region with spontaneous Ising order), and 
hence the overall amplitude of (5.2) is determined only 
by (4vcn) + (rfodd)- We show our numerical results for 
these, and o ther related quantities, for the Hamiltonian 
TL d in fll.l2j ) in Figure y. There is a broad maximum 
in (devon) + (d dd) near the Ising critical point, as this 
is the region with maximal dipole number fluctuations. 
The critical singularity in this quantity at A = A c is de- 
termined by that of the energy operator of the Ising field 
theory: this singularity is weak and is essentially unob- 
servable in Fig [l4|. The quantities sensitive to the Ising 
order parameter (such as (d even ) — (e? dd)), show more 
singular behavior in Fig [l4| near A c determined by the 
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FIG. 14: Ground state exp ectat ion values of {de) and {d\di) 
for the D — 1 model TC d in ( |i.l2| ). The results are for TV = 16 
sites and periodic boundary conditions. A very small ordering 
field was applied to choose one of the degenerate Ising ground 
states present for sufficiently negative A. We have chosen the 
gauge in which (d) are real. 



N ^ y ' 



dldi 



(5.3) 



in D = 1 (sec (2J3), and a related definition can be made 
for higher D. One possibility for coupling to the Ising 
order parameter experimentally would be to introduce 
a phase-locked subharmonic standing wave at half the 
wave vector of the optical lattice so that the standing 
wave takes the form (in ID) 

$(x, t) oc [cos(Qx) cos(Qct) + B cos(Qa;/2 + 0) cos(Qct/2)]. 

(5.4) 

Squaring this and taking the time average gives the ef- 
fective lattice potential 



V(x) oc — [cos {Qx)+B 2 cos 2 (Qa;/2 + i 



(5.5) 



Adjusting the relative phase to 9 = or tt/2 adds a 
'staggered magnetic field' term to the Ising Hamiltonian 



H B oc ±B 2 



(5.6) 



magnetization exponent /?. However these observables 
are not detectable by a measurement of the momentum 
distribution function. 

In higher dimensions (D > 1) for the case where e is 
aligned along one of the lattice directions, the dependence 
of the distribution function on qu should be qualitatively 
similar to the q dependence in the D = 1 case discussed 
above. However a much clearer signal of the transverse 
superfluidity should be visible. The presence o f th e (p) 
and (h) condensates imply that the correlator ( |5.l| ) has 
phase coherent contributions when rj— r& lies in the plane 
perpendicular to the applied 'electric' field. This implies 
that in states with transverse superfluidity, there should 
be Bragg peaks along lines in q space with values of 
equal to the reciprocal lattice vectors of the D — 1 dimen- 
sional lattice lying in the plane perpendicular to e. As the 
transverse dimensionality is D — 1 = 2, the superfluid or- 
der can only be quasi-long-range at nonzero temperature, 
and hence the Bragg peaks are not true delta functions 
in the infinite volume limit, but are power-law singular- 
ities. Experimental detection of these Bragg lines would 
be quite interesting. 



A simpler experimental method for the case where the 
trap confinement is strong in the directions transverse to 
the axis of the ID lattice is the following. An additional 
standing wave (derived from the same laser) but oriented 
in the y direction (say) would yield 

$(x, y = 0,t)(x [cos(Qx) cos(Qct) + B cos(Qct)} (5.7) 
and hence a potential along the x axis of 

V(x,y = 0) oc -[cos 2 {Qx) + 2B cos(Qx) + B 2 } (5.8) 

which would also couple to the Ising order 

H B oc Bcf>. (5.9) 

In either case, such a perturbation could be used to break 
the Ising symmetry and selectively populate one of the 
two Ising states. In addition, it could be used to mea- 
sure the order parameter itself. The AC stark shift of 
the atomic hyperfine levels would differ between adja- 
cent sites. The relative strengths of the split hyperfine 
absorption lines .Mould then be a measure of the Ising 
order parameter £3 



B. Ising Order Parameter 

We have seen that the Ising order is not directly re- 
flected in the momentum distribution and hence can 
not be measured in the free expansion method described 
above. The properties of the Ising order parameter, 
are discussed in Appendix H; one convenient definition 
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APPENDIX A: FLUCTUATIONS AND 
QUANTUM FIELD THEORIES: 
SUPERFLUID-INSULATOR TRANSITIONS 



layer. The values of \xe are to be determined at the end 
by the requirements 



The mean field theory of Section [II B can be used 
as a starting point for a more sophisticated treatment 
of fluctuations. Such fluctuations will modify the mean- 
field exponents in the vicinity of the second order phase 
boundaries in Fig|Io|a, but could also change the topology 
of the phase diagram to that in Fig [To|b. 

We analyze fluctuations about the mean field results 
using a method very similar to that described in Chapter 
10 of Ref. [ll] for the Hubbard model. We decouple the 
intraAayer hopping terms in TL p h (those in the last line 
of ( 1. 18| ) oruv by Hubbard-Stratonovich transformations 
using complex fields Pe (r j_ , t) and He (r j_ , r) where is 
a spatial co-ordinate for the D — 1 transverse directions, 
and r is imaginary time. Then, after standard simplifica- 
tions, we obtain an expression for the partition function 
Z p h of Tlph which has the following schematic form 



dhxZ 



pit 



Z ph = J VP t (T U T)VH t (T X ,T) 

exp - J d z? - 1 r ± (5 +5 1 ) 



(Al) 



The action iSo involves couplings only within a single 
layer £, but with different values of rj_ 

So = [ drJ2 [K P \V x P t {v U T)\ 2 +t p \P 1 {v u t)\ 2 

+K h \V ± H e (r ± ,T)\ 2 +r h \H i (r ± ,T)\ 2 ] , (A2) 
and K p h, r p h are coupling constants. Note that the 



factors of no and no + 1 in the last line of ( l.lq ) break 
particle-hole symmetry and so there is no special sym- 
metry relation between these coupling constants. The 
action Si couples different layers and times together for 
the same value of 



V Pt (T)Vhe(T)P\pe(T),h e (T)] 



x exp 



- / dr 



+TL P h, 




with TL p h,mf defined in (3.7), and V is a pr ojection oper- 
ator which represents the constraints (3.8) (these could 
be imposed formally in the functional integral by a very 
strong on-site repulsive inter action among the pe and hg 
bosons) . As i n Section III B , we have imposed the con- 
straints (1.17) by time-independent Lagrange multipliers 
("chemical potentials") fig; as we noted earlier, there is 
no approximation involved in neglecting the fluctuations 
of fit, because there is only one constraint per layer and 
there are a macroscopic number of particles within each 



= 



(A4) 



Further progress in describing the properties of Z p h re- 
quires some understanding of the structure of <Si. This 



was already addressed to some extent in Section III B 
where we explored the properties of the Hamiltonian 
TL p h,mf- However, here we need to generalize that anal- 
ysis to the case where its arguments are time-dependent 
fields Pe(r±, r) Hi(r±,r). This is quite an involved task, 
but we will only need some general constraints that are 
placed on the structure of Si by the principles of gauge 
invariance. In particular, associated with the conserva- 
tion laws (1.17), we observe that Z p h is invariant under 
the time and layer-dependent transformations generated 
by the arbitrary field 4>e(j) 



Pe+i 




hi 


-► /i £ e~^ (T) 


Pt+i 


-> P e+ ie lMT) 


H e 






. -®4>e 


fie 


-> fie + i^- 

OT 



(A5) 



We are interested here only in the case of time- 
independent and so this transformation takes /j,£ into 
an unph ysic al set of values; nevertheless, as we will see 
shortly, (A5) is still useful in placing constraints on Si in 
the physical regime. 

First, we address the influence of fluctuations by ap- 
proaching the transition involving condensation of Pe, He 
from the side of large and positive A, i.e. we increase E 
(and decrease A) until mean-field theory indicates we are 
approaching a phase with transverse superfluidity at the 
point A in Fig |l^. The ground state of TL p h, m f is trans- 
lationally invariant in this region, and so we can safely 
assume that all the coupling constants in Si are also in- 
dependent of £. Similarly, we can assume that m is also 
independent of I. If we were to approach the condensa- 
tion of Pe, He from the opposite side of negative A, the 
ground state of H. p h. m f would have a broken Ising sym- 
metry, and the following analysis would only need to be 
modified by allowing all couplings, and /i, to depend upon 
the £ sublattice. We describe the action Si by expanding 
it in powers of the fields Pe, He, and in their temporal 
gradients (the rj^ and r dependence of these fields is now 
implicit); to second order in the fields and to first order 
in temporal gradients, the most general terms invariant 
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under (A5) are 



Si 



,dP £ 



.dH, 



K ph I Pt 



dH, 



e+i- 



(A6) 



8t 

+ r p \P e \ 2 + 7 h \H e \ 2 + r ph (P e+1 H e + P* l+1 H^) 



Consistently requiring invariance of (]A6|) under the time- 
dependent gauge transformations ( |A5| ) to the order we 
have performed the expansion in S\ demands additional 
constraints on the coupling constants above; these are 



K 

p ~ dp 



K h = 



dr h 
djjb 



K p h — 



dr ph 
9/i 



(A7) 



There are also a large number of permitted higher order 
terms in S\ which we have not written down explicitly; 
some of these will play an important role below. 

Armed with the low order terms in the action Sq + Si 
controlling the fluctuations of Pg and Hi we can now 
use standard techniques to focus on the low energy exci- 
tations. It is natural to diagonalize the quadratic form 
displayed in these actions: this will lead to two eigen- 
modes with distinct eigenvalues. We focus attention on 
the lower eigenmode, while integrating out the higher 
eigenmode. We identify the lower eigenmode by the field 
*f. this has the structure 



Vi(r ± ,T) = CpPi +1 (T_L,T) + C h Hg{T±,T) 



(A8) 



for some constants Cp t h- Note that we are performing 
the same 'rotation' in field space for all and r (and 
hence all frequencies) . This ensures that *i has a simple 
behavior under (|A5r): 



i<t>t{j) 



(A9) 



We integrate out the high energy mode orthogonal to 
( |A8| ) , and obtain our final effective action now expressed 
in terms of *f. 



Cn = IV , ^ 



r>|^| 2 + ~|^ 



+ v\* e \ 2 \* e+1 \ 2 



dr 



(A10) 



We have rescaled *e and time to obtain unit coefficients 
for the first two terms in £ > and the rj_ dependence of 
*i is implicit. We have also written down a quartic non- 
linearity within a layer (it), and the simplest coupling 
between neighboring layers (v) which preserves invari- 
ance under (A9); we expect both these couplings to be 



positive because of the repulsive interactions between the 
microscopic bosonic degrees of freedom. The parameter 
Tip tunes the system across the quantum phase transition 
at the point A in Fig [w] which resides at = r^ c \ the 
transition is from the featureless, gapped phase at large 
positive A (r^ > r^ c ) to a phase with superfluidity in the 
transverse D—l dimensions as A is decreased (r^, < r^, c ); 
the superfluidity is associated with the condensation of 
* f . 



Just as in the derivation of (A7), we can also exam- 
ine the cons equence s of t ime- dependent gauge transfor- 
mations in (A5.A9) on (A1C). This now leads to the 
relationship 



(All) 



Combined with (A4), the above result now yields a cru- 
cial result. Close to the quantum critical point, the singu- 
lar free energy associated with Z p h is d eter mined directly 
by r,/,. For this singular term to obey (A4), we conclude 
(as also argue d in C hapter 10 of Rcf[ll]) that dr^/dji = 



at r,h 



(All) now implies 



— atr^i = r^ c 



(A12) 



So we can neglect £ 1; and the critical theory is described 
entirely by Cq. Within each layer I, this theory has the 
relativistic invariance of [D— 1) + 1 spacetime dimensions 
and dynamic critical exponent z = 1. 

Be fore turning to an examination of the properties of 
(AlC), we pause to discuss the modifications required 
to describe the onset of transverse superfluidity with in- 
creasing A in the region A < at the point B in Fig 10 (a 
similar reasoning can also be applied to the point C2 in 
Fig [lO^) ■ Here, long-range Ising order is already present 
in H p h m f for A sufficiently negative. We can proceed 
to a description of the superfluid transition as above, 
but as noted earlier, all couplings in (A10) will acquire 
an £ dependence which modulates with period 2. The 
tuning parameter will also be different for even and 
odd £. Consequently only *£ with I even (say) will be- 
come critical near the transition, while *i with I odd 
remains non-critical and can be integrated out. The sim- 
plest interlayer coupling between critical modes is now 
■i+2\ 

likely to be attractive. 



We now return to an examination of ( A10) for the case 
of ^-independent couplings at the transition with U — E 
positive at the point A in Fig It remains to exam- 
ine the consequences of the interlayer coupling v on the 
standard theory of the superfluid-insulator transition. At 
v = 0, we have the standard (p 4 field theory with 0(2) 
symmetry in (D — 1) + 1 = D spacetime dimensions. As 
a first step, we can compute the scaling dimension of v at 
its critical point. A standard power-counting argument 
shows that 



2 a 

dimful = D = — 

v v 



(A13) 
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where v and a are the standard correlation length and 
"specific heat" exponents in D spacetime dimensions. In 
D = 3, the 0(2) fixed point hasta a = -0.015 < 0, and 
so we conclude that v is formally irrelevant. In I? = 2, 
the very weak specific heat singularity at the Kosterlitz 
Thouless transition suggests the same conclusion. 

A more complete analysis of the influence of v can be 
obtained by considering a physical susceptiblity for or- 
dering in the lo ngitu dinal direction. As we have seen in 
Sections [II] and III B , the simplest allowed ordering is a 
density wave of period two. The tendency to this order- 
ing is measured by the static susceptibility x-n 



-T 

11 £,£' 



<|vMr ± ,r)| 2 |vMO,0)| 2 ). (A14) 



Note that this response function is similar to S n in fl2.3| ) , 
but we are considering here a zero frequency response, 
while ( 2.3) inv olved an equal time correlator. We can 
compute (A14) in powers of v, and by a familiar Dyson- 
type argument, write it as 



X* = 



1 - 2vC 



(A15) 



where C is an 'irreducible' correlator within a single layer 
(it is irreducible with respect to cutting a v interaction 
line): 



C 



jD-1 



r ± dT(\y e (r ± ,T)\ 2 \V e (0,0)\ 2 ) . (A16) 



The computation leading to ( A15] , A16) is the field- 
theoretic analog of the computations which lead to a 
dipole bound state induced by the inte rlayer coupling 
in the strong-coupling analysis of Section III A . Ignoring 
the influence of v on C, standard scaling arguments imply 
that C has a singular part which behaves as 



C 



(A17) 



If we had a > 0, then the denominator in (A15) would 
vanish at some > r^ c for any small v, and x* then di- 
verges: this would imply the presence of an Ising density 
wave transition before the onset of superfluidity. How- 
ever a < in D = 3, and so this condition does not 
apply. Nevertheless, there is a significant (albeit finite) 
enhancement of the specific near the 0(2) critical point 
in D = 3, and so the instability in x-k may well occur 
for a moderate value of v. If so, the mean-field phase 
diagram would be modified, and the Ising ordered phase 
would fully overlap and extend beyond the region with 
transverse superfluidity. Indeed, under suitable condi- 
tions, the superfluid phase could also shrink to zero, and 
we would then have only a single Ising transition be- 
tween two insulating phases. Alternatively, if the Ising 
fluctuations are weaker, Xn could diverge somewhere in 
the superfluid phase to the left of A in Fig |l^, and then 
the mean-field phase diagram would be modified to the 
structure in Fig UX 



APPENDIX B: ISING PHASE TRANSITION 

In Appendix ^ we completed a description of fluc- 
tuations near all the superfluid-insulator transitions in 
Fig [l(| It remains to describe the second order Ising 
critical point C\ in Fig [lO|b; this we do in the present 
appendix. 

As usual, we expect the Ising phase transition to be 
realized by a quantum field theory of a real scalar field 
(f>(r, t), where r = (^, r jJ is a D dimensional spatial 
co-ordinate. The main subtlety here is that the Ising 
transition occurs in a background of transverse super- 
fluid order, and corrections from superflow fluctuations 
can lead to anisotropic singular corrections to the criti- 
cal theory. A theory of an Ising order parameter coupled 
to isotropic superfhpw fluctuations has been analyzed by 
Frey and BalentsLJ; here, we will show that the partic- 
ular anisotropic nature of both the superfluid and Ising 
order leads to a more singular coupling between the two 
order parameters. 

Any observable sensitive to the period 2 modulation 
in the density of particles or holes can be used to define 
the order parameter <fi(r, r). A convenient choice in our 
present continuum formulation is to take 



,r x ,T)~(-l)^(rL,T)| s 



(Bl) 



An effective action, S$, for the Ising field 4> can be gen- 
erated by using as a Hubbard-Stratonovich field to de- 
couple the v term in ( A1C| ) . This leads to an action with 
the structure 



d D rdr 



+ ^(V||0) 2 +^ 



U'l 



J2 I d D -VdT(-l)^(r ± ,T)| 2 ^,r ± ,T), 

(B2) 



where the fluctuations of are described by (A1C), and 
we have included the usual analytic terms present in the 
4> 4 theor y of an Ising quantum critical point. The last 
term in (B2) represents a linear coupling between the 
Ising order parameter and density fluctuations in the su- 
perfluid state. In the isotropic case considered by Frey 
and Balents such a linear coupling was absent, and the 
simplest allowed coupling was between <p 2 and the density 
fluctuations: this was because the Ising order parameter 
represented a density wave at a large wavevector, and 
they coupled linearly only to fluctuations of the super- 
fluid phase at the same wavevector, and the latter are 
quite high energy. In the present case also, the (— 1) 
factor in the last term in (p2[), also shows that <f> couples 
linearly to the superfluid phase fluctuations at a wavevec- 
tor q\\ = 7T. However, the key difference here is that the 
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superfluidity is present only along the transverse direc- 
tion, and, to leading order, the superfluid phase fluctua- 
tions are independent of qu 



The singular effect of the wi term in (B2) can be illus - 
trated by integrating out the using the action (AlC) 
in a single loop approximation. To leading order in u, we 
are in the transverse superfluid state as long as < 0, 
and a simple calculation of the phase and amplitude fluc- 
tuations of the superfluid order parameter shows that we 
generate the following term in <S^: 



9ii ,qj-,w 



(B3) 



where u> is an imaginary frequency. Note that this is a 
singular function of qx and uj only when ^ 0. We 
do not expect K$ = near the Ising critical point, be- 
cause exact particle-hole symmetry is not present in the 



un derly ing Hamiltonian, and the arguments which lead 
to (|AlSf ) hold only at the superfluid-insulator transition. 
All the couplings in 5^ can be expected to be a smooth 
function of [i, and the constraint-is now expected to lead 
only to a Fisher renormalizationE3 of exponents. An anal- 
ysis of <S^ with (B3) included requires a renormalization 
group computation: this we leave to future work, as a 
full discussion of the renormalization of the momentum 
dependence of the propagator requires a two loop analy- 
sis. 

In closing, we note that although ^ 0, in practice 
the degree of particle-hole symmetry breaking is quite 
small, as is indicated by the almost equal values of (p) 
and (h) in Figs |l| and ||. So can also be expected 
to be quite small, and we shoul d, th erefore, also consider 



the case = 0. In this case, (B3) docs not induce any 



singular terms, and we have to consider terms induced by 
fluctuations at higher orders in u, and also the term 
included in Ref [L9l 
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